load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"


;#HERE YOU WILL PUT THE PATH 
location="/Users/mypath"

filsw= systemfunc ("ls "+location+"lowice/siconc_allmon_ymonmean*.nc")
filsc= systemfunc ("ls "+location+"control/siconc_allmon_ymonmean*.nc")
f_wgt = addfile(location+"/mesh_mask_nemo.Ec3.2_O1L75.nc","r")
fw    = addfiles (filsw, "r")
fc    = addfiles (filsc, "r")

ListSetType (fw, "join")
ListSetType (fc, "join")
print(fw)
print(fc)

;*******************************************;
;Read in variables
;*******************************************
lat2d=fc[0]->latitude
lon2d=fc[0]->longitude
nlat=292
nlon=362
printVarSummary(lat2d)
 
   ens=10  

   londis=f_wgt->e1t(0,:,:)
   printVarSummary(londis)
   latdis=f_wgt->e2t(0,:,:)
   printVarSummary(latdis)
  
    wgt=londis
    wgt=londis*latdis

   steps=12  ; 12 months

tempcss =  new((/ens*steps,292,362/),float)
icesumc =  new((/ens,steps/),float)

do i=0,ens-1

time=fc[i] ->time
ntime= dimsizes(time)

tempsc = fc[i]->siconc

do k=0,291
do j=0,361

if (lat2d(k,j).lt.0) then
tempsc(:,k,j)=0
end if
end do
end do


icesumc(i,:) = wgt_areasum2(tempsc(:,:,:),wgt(:,:),0)
icesumc(i,:) = icesumc(i,:)*10^(-14)

printVarSummary(icesumc)
delete(time)
delete(ntime)
delete(tempsc)
end do

print(icesumc)


icesumw =  new((/ens,steps/),float)

do i=0,ens-1

time=fw[i] ->time
ntime= dimsizes(time)

tempsw = fw[i]->siconc

printVarSummary(tempsw)

do k=0,291
do j=0,361

if (lat2d(k,j).lt.0) then
tempsw(:,k,j)=0
end if
end do
end do


icesumw(i,:) = wgt_areasum2(tempsw(:,:,:),wgt(:,:),0)
icesumw(i,:) = icesumw(i,:)*10^(-14)

printVarSummary(icesumw)

delete(time)
delete(ntime)
delete(tempsw)

end do

;***********************************
; plot
;***********************************
wks = gsn_open_wks("ps","icearea_monthly_ecearth")


  month = ispan(0,11,1)*1.0
   plot    = new (2,"graphic")


   res                      = True          ; individual plot
   res@gsnDraw              = False
   res@gsnFrame             = False
   res@xyLineThicknessF     = 2.0
   res@xyMonoLineColor         =True
   res@xyLineColor         = "black"

   res@tmXBMode        = "Explicit"
   res@tmXBValues      = month
 
   res@tmXBLabels      = (/"J","F","M","A","M","J"\
              ,"J","A","S","O","N","D"/)
   res@xyDashPattern        = 0                  ; Make curves all solid

   res@vpHeightF= 0.4                    ; change aspect ratio of plot
   res@vpWidthF = 0.8

   res@trXMinF  =  0
   res@trXMaxF  = 11
   res@trYMinF=0
   res@trYMaxF=16

   res@tiXAxisFontHeightF = 0.010
   res@tiYAxisFontHeightF = 0.015

   res@gsnLeftString = ""
   res@tiXAxisString = ""
   res@tiYAxisString = ""

   icesumc!0 = "ncl_join"
   icesumc!1= "time"
   icesumc&time= month

   icesumw!0 = "ncl_join"
   icesumw!1= "time"
   icesumw&time= month


   plot(0)  = gsn_csm_xy (wks,month,icesumc(0,:),res)

   do i=0,8
   plot0  = gsn_csm_xy (wks,month,icesumc(i+1,:),res)
   overlay(plot(0), plot0)
   delete(plot0)
   end do

   res@xyLineColor         = "cyan3"
   do i=0,9
   plot1  = gsn_csm_xy (wks,month,icesumw(i,:),res)
   overlay(plot(0), plot1)
   delete(plot1)
   end do

 ;***** make panel *****
   resP                     = True          ; panel resources
   resP@gsnMaximize         = True

   resP@txFontHeightF = 0.015
   gsn_panel(wks,plot,(/1,1/),resP)

